###########################################################################
####State Capacity, Insurgency, and Civil War - Grid Cell RE Robustness####
###########################################################################
library(MASS) 
library(pscl)
library(foreign)
library(Hmisc)
library(memisc)
library(Matrix)
library(lme4)
library(lattice)
library(gamlss.mx)
library(doBy)
library(car)
library(plotrix)
library(glmmADMB)
library(MCMCglmm)
library(stargazer)
library(mvtnorm)
library(Zelig)
library(readstata13)
library(pROC)
library(robust)
library(DAAG)


########################
###Create The Dataset###
########################
###Importing the main data
#set working directory to location of dta file
setwd("~/Google Drive/NSLC/Final Dofiles/")
#Import the data
main.data <- read.dta("ntlgrid.dta")

###Model 1: Civil War
#subset the data
dat.1<-na.omit(main.data[,c("onset", "lnNL_sum", "lagcivconflagtemp", "loglagpop", "loglagppp", "mnt1", "logbdist1", "logross_oil_prod", "nwstate", "instability", "lagp_polity2", "ethfrac", "relfrac", "logcellarea", "logttime", "year93", "year94", "year95","year96", "year97", "year98", "year99", "year00", "year01", "year02", "year03", "year04", "year05", "year06", "year07", "year08", "gid")])

#Create prior
zi.prior <-  list(R = list(V = 1, n = 1.002, fix=1),
                  G = list(G1 = list(V = 1, n = 0.002)))

#Model 1: civil war
MRE.1.o.c <- MCMCglmm(onset ~ lnNL_sum + lagcivconflagtemp + loglagpop + loglagppp + mnt1 + logbdist1 +  
                        logross_oil_prod + nwstate + instability + lagp_polity2 + ethfrac + relfrac +
                        logcellarea + logttime + 
                        year93 + year94 + year95 + year96 + year97 +
                        year98 + year99 + year00 + year01 + year02 + year03 + year04 +year05 +
                        year06 + year07 + year08, 
                      data = dat.1, family = "categorical",
                      random =  ~ gid,
                      prior = zi.prior, slice = T,
                      nitt = 210000, burnin = 200000, thin = 10,
                      verbose = TRUE, pr = TRUE, pl = TRUE)
summary(MRE.1.o.c)

###Model 2: Ethnic War
#subset the required data for ethnic war
main.data.eth <- main.data[ which(main.data$second > 0.04999),]
#subset data for analysis
dat.2<-na.omit(main.data.eth[,c("onset", "lnNL_sum", "lagcivconflagtemp", "loglagpop", "loglagppp", "mnt1", "logbdist1", "logross_oil_prod", "nwstate", "instability", "lagp_polity2", "ethfrac", "relfrac", "logcellarea", "logttime", "year93", "year94", "year95","year96", "year97", "year98", "year99", "year00", "year01", "year02", "year03", "year04", "year05", "year06", "year07", "year08", "gid")])

#Create prior
zi.prior <-  list(R = list(V = 1, n = 1.001, fix=1),
                  G = list(G1 = list(V = 1, n = 0.001)))

#Model 2: ethnic war
MRE.2.o.c <- MCMCglmm(onset ~ lnNL_sum + lagcivconflagtemp + loglagpop + loglagppp + mnt1 + logbdist1 +  
                        logross_oil_prod + nwstate + instability + lagp_polity2 + ethfrac + relfrac +
                        logcellarea + logttime + 
                        year93 + year94 + year95 + year96 + year97 +
                        year98 + year99 + year00 + year01 + year02 + year03 + year04 +year05 +
                        year06 + year07 + year08, 
                      data = dat.2, family = "categorical",
                      random =  ~ gid,
                      prior = zi.prior, slice = T,
                      nitt = 210000, burnin = 200000, thin = 10,
                      verbose = TRUE, pr = TRUE, pl = TRUE)
summary(MRE.2.o.c)

###Model 3: Civil War
#subset the data
dat.3<-na.omit(main.data[,c("onset", "lnNL_sum", "lagcivconflagtemp", "loglagpop", "loglagppp", "mnt1", "logbdist1", "logross_oil_prod", "nwstate", "instability", "laganocracy", "lagpolbin", "ethfrac", "relfrac", "logcellarea", "logttime", "year93", "year94", "year95","year96", "year97", "year98", "year99", "year00", "year01", "year02", "year03", "year04", "year05", "year06", "year07", "year08", "gid")])

#Create prior
zi.prior <-  list(R = list(V = 1, n = 1.002, fix=1),
                  G = list(G1 = list(V = 1, n = 0.002)))

#Model 3: civil war (binary democ)
MRE.3.o.c <- MCMCglmm(onset ~ lnNL_sum + lagcivconflagtemp + loglagppp + loglagpop + mnt1 + logbdist1 +  
                        logross_oil_prod + nwstate + instability + ethfrac + relfrac +
                        laganocracy + lagpolbin + logcellarea + logttime +
                        year93 + year94 + year95 + year96 + year97 +
                        year98 + year99 + year00 + year01 + year02 + year03 + year04 +year05 +
                        year06 + year07, 
                      data = dat.3, family = "categorical",
                      random =  ~ gid,
                      prior = zi.prior, slice = T,
                      nitt = 210000, burnin = 200000, thin = 10,
                      verbose = TRUE, pr = TRUE, pl = TRUE)
summary(MRE.3.o.c)